Phase diagram of 3D ANNNI model in an effective-field approximation 
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An efTective-field method for calculation of thermodynamic properties of three-dimensional lattice 
spin models is developed. It is applied to the ANNNI model on the simple cubic lattice. The phase 
diagram of the model, consisting of a large number commensurate phases and of an incommensurate 
phase, is calculated, confirming the results of previous approaches. The phase transition lines for 
a number of commensurate structures are localized and a strong evidence for absence of the direct 
phase transition between commensurate phases and the disordered phase is found. 

PACS numbers: 64.70.Rh 



I. INTRODUCTION 



In this paper we study the ANNNI model on a simple 
cubic or tetragonal lattice. This model was first intro- 
duced by Elliott [H in order to understand modulated 
magnetic materials. It is reviewed by Selke and Yeomans 

Ij, . The model is known to form a low temperature 
ferromagnetic phase for a small next-nearest-neighbour 
interaction and a (2) phase for a large one. The wedge in 
the nnn interaction-temperature phase diagram between 
this two phases is, at low temperatures, filled by infinite 
number commensurate phases. 

Theoretical study of the ANNNI model has been 
based on a large number of various approaches. The 
devil's staircase structure of the phase diagram at 
low and medium temperatures was elucidated by low- 
temperature series expansion |^ and mean-field approx- 
imations 0) Q) S 0) US • The incommensurate phase was 
also treated by the free fermion approximation pH IT^ . 
Recently, an anisotropic scaling at the Lifshitz point was 
used to calculate several critical exponents at this point 
[l^. A considerable effort was also devoted to investiga- 
tion of ANNNI thin films Ellilll^. 

The mean-field approximations describes qualitatively 
well the phase diagram of the ANNNI model, neverthe- 
less, some of its features were challenged by other ap- 
proaches, e.g. the stability of the commensurate phase 
up to the transition line to the disordered phase. 

To improve the performance of mean-field treatment of 
the ANNNI model, we develop an effective field method, 
which is a generalization of the cluster transfer-matrix 
method successfully applied to 2D spatially modulated 
structures pl.lTlll6|. 

Our effective field method resembles the nonlinear 
mapping approach of Bak 0,112, but, instead of magne- 
tization, it maps a large number of effective fields, which 
simulate the cluster environment. It is related also to the 
DMRG method [H, and for the 2D ANNNI model they 
yield similar results [isl IT^ . Comparing with DMRG 
approach, our method is much simpler, and instead of 
diagonalization of density matrix and renormalization of 
transfer matrix by matrix multiplication it requires only 
calculation of square-root of a function of cluster spin 
configurations and real- number multiplications 23] . The 
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FIG. 1: Spin-spin interactions in 3D ANNNI model; nn in- 
teractions Jo and Ji are ferromagnetic, nnn interactions J2 
are antiferromagnetic. 



results of our method is in general agreement with other 
approaches, and it removes the artefacts of the previous 
mean-field methods. 

In Section II the 3D ANNNI model is specified, and a 
new effective field approximation is developed. Results 
of numerical calculations and a tool for distingushing be- 
tween commensurate and incommensurate phases, which 
lead to construction of the phase diagram are presented 
in Section III. 



II. MODEL AND METHOD 

We shall generalize the cluster transfer-matrix method 
(an effective-field approximation) developed and applied 
to 2D space- modulated structures some time ago ^3 
.16, 23J. The development of the 3D method follows the 
same ideas that were used in 2D case, however, a new 
approximation in the course of calculation has to be done. 
For reasons of clarity the method is developed for an 
ANNNI-type model but it can be easily reformulated for 
any 3D model with short-range interactions. 

The three dimensional ANNNI model on a simple cubic 
lattice consists of two dimensional planes, within which 
each spin is coupled to its nearest neighbors by a ferro- 
magnetic interaction Jq. However, in the direction per- 
pendicular to the planes, spins are coupled by compet- 
ing ferromagnetic nearest- and antiferromagnetic next- 
nearest-neighbor interactions (Fig. 1). For reasons of 
simplicity Ji = Jq is further assumed. 
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As the interactions between spins cr^ j fe 
Hamiltonian of the 3D ANNNI model 



±1 in the 



H = ^ ~JQcri,j,k{cri+i,j.k + o'ij-+i,fc) — Ji(yi,j,kcrij,k+i 

+ J2'^i,j,kO'iJ,k+2 (1) 

involve only three layers, it can be written as a sum of 
layer Hamiltonians Hi(Si, Si+i, Si+2) which depend on 
three layer variables Sk = {uij^k} hj G (—00,00). Since 
there are only nearest-neighbor interactions inside the 
layers, the layer Hamiltonian can be expressed as a sum 
of cluster Hamiltonians defined on 2 x 2 x 3 clusters with 
the longer side oriented along the J2 interaction. 

H = '^^{Hk{Sk, Sk+l, Sk+2) — ^ ^ Hij^kio'l,m,n) 
k k i,j 

(2) 

where Z = i, z + 1, m = j,j + l, n = k,k + l,k + 2. 

The exponential of the layer Hamiltonian 
Hk is further denoted by Tk{Sk, Sk+i, Sk+2) = 
exp[(3Hk{Sk, Sk+i, Sk+2)] and sometimes called transfer 
matrix though it is rather a function of spin variables. 

Then the summation in the partition function 

Z=J2 exp[/3i^(aO] 

may be performed consecutively layer by layer generating 
a set of auxiliary functions ^fc and normalization factors 

kiSk, Sk+l)Tk{Sk, Sk+1, Sk+2) — ^k'^ k+l{Sk+l, Sk+2) 

(3) 

starting from an appropriate function ^'i(5i,S'2) that 
may be interpreted as a boundary condition of the system 
on a semi-infinite lattice. The values of ^'^ for fc — > 00 
mostly do not depend on the input except the vicin- 
ity of a first order phase transition. Here the different 
bulk values correspond to one stable and one or more 
physically unstable solutions. The stable solution is the 
one with the lowest free energy that is proportional to 

As we see, the auxiliary functions in the transfer ma- 
trix method are some general positive functions defined 
on clusters of planes in 3D models. For lower-dimensional 
models they are defined on clusters of rows in 2D and 
clusters of sites in ID. In a one-dimensional model, the 
auxiliary functions depend on finite number of spin vari- 
ables, in 2D and 3D cases they acquire infinite number 
of values which cannot be generally found by numerical 
calculations. Instead of the whole function at the right- 
hand side of (3), we further calculate only its correlation 
function, sum of '^k+iiSk+i, Sk+2) over the whole lattice 
except a small cluster, and the true auxiliary function ^' 
is approximated by a more convenient one, nevertheless, 
exactly reproducing the calculated correlation functions. 

As we do not use any further information from the left- 
hand side of Eq. (3), all the remaining properties of the 



approximate function are derived from the requirement 
of maximum of the information entropy S* = ^ log 21] . 
To maximize S under the condition that the partial sum 
of ^ is equal to the given correlation function, Lagrange 
multipliers corresponding to each configuration of the 
cluster have to be introduced. It is easy to show that the 
desired auxiliary function can be expressed as a product 
of exponentials of the cluster Lagrange multipliers. Thus, 
the requirement of maximum of the information entropy 
leads to a factorization of the auxiliary function. If only 
factorized auxiliary functions are used, the left-hand side 
of (3) is completely factorized for short-range interactions 
and its partial summation is equivalent to calculation of 
a correlation function of a statistical system of the di- 
mension lower by 1 than that of the original problem. 
It means that for 2D system this step can be performed 
exactly, but for 3D the factorization procedure must be 
applied even for calculation of the correlation function. 

In the case of 2D model (ID auxiliary functions) the 
application of the above considerations is straightfor- 
ward. Let us denote the correlation function on a small 
cluster of the length n + 1 by 



i-\-7i 



- E 

3e(-oo.i)(i+7, 



^i'^j), 



which is assumingly known from the previous calculation 
step. The approximate function ^ is defined as a product 
of unknown cluster functions 0^^" {aj ) defined on n -|- 1 
sites: 



^^l[e^+"{a,) Je{^,^ + n}. 



(4) 



We would like to express the cluster functions 0*^" in 
terms of $^+". 

Let us denote the left eigenfunction of the function 
(transfer matrix) 0*^" 



Y.ei^-\a,)e^-{a,)^X6^{a,) 



(5) 



by 6 and its eigenvalue by A. 



9j+"-i and 91+"^ are iden- 
tical function defined on different clusters if we do not 
expect any space modulation in this direction.) 

Since is the result of summation of ^ from —00 

to i, correlation function corresponding to ^ can 

be expressed as 



i-\-n 



2+n — l0?+n^?+l 



where 9 is the right eigenfunction of B defined by 

Since we require = the unknown cluster func- 

tion 6!+" is 



Q^+n ^ 



i+n 



ri+n-lni+1 ■ 
'i "i+n 



(7) 
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Unfortunately, the eigenfunctions and 9 are implicit 
functions of O. On the other hand, it can be easily shown 
that 



ni+n— 1 



e 



■ e 



i+n— 1 




(8) 



have the same eigenvalues as the original cluster functions 
e^+". Substituting (7) for 9,^+" we get 
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ri+n—lni ni+nni+l 



where 



(9) 

= E.. <fr" = ^S'l ^Itn and similarly 
^. Thus, obeying the condition of maximum of 
information entropy, relations (4, 9) yield a possibility 
to express the approximate chain auxiliary function ^ in 
terms of the known correlation function 



i-\-n 



i i+1 

(10) 

In the case of 2D auxiliary function (for 3D models) 
the relation (10) is further valid, only the indices denote 
infinite rows of sites rather than sites. However now, we 
cannot expect that the correlation function <!)*"'"" defined 
on n + 1 infinite rows could be found in previous cal- 
culations, but only a function on a finite cluster of the 
size n X I. We denote it by where the first in- 

dices represent rows and the second ones columns of the 
lattice. The most serious difference between 3D and 2D 
models is that the eigenfunction in (5) cannot be found 
exactly but it must be factorized, as well. Then all the 
functions in the expression for <!)*"'"" are factorized simi- 
larly as at the left-hand side of Eq. (3), i.e. the whole 
procedure that we applied to the chain auxiliary function 
^I^, and that has led to (9), can be applied to n-row and 



i+n — l 



n - 1-row correlation function and $^ 

respectively, appearing in 2D version of (10). We obtain 



n 



i~\-n,j+l 



etc. 



(11) 



Thus, by consecutive application of the factorizing pro- 
cedure (11) to all terms in (10), an approximation \1/ to 
the function 5* can be expressed in terms of its cluster 
correlation function 



' hi 



(12) 



fce(- 



=o,i)U(i + »i,oo) 
oo,3)U(j + i,oo) 



The expression reads 

* « * = a., 



(13) 



n 



^i + n.3 + i 4/^> + n-l,j + i-l^, + „-l.j + ;^, + „,3 + i-l^i + ,i,j + l 



^,3 



Unlike in ID case, the correlation function calculated 
from ^' is only approximately equal to that calculated 
from It would be true if we were able to factorize 
the whole two-dimensional plane function O^'"*"" in 2D 
version of (10) and not only each correlation function $ 
separately. 

All the functions in (13) are plane dependent in the 
case of a modulated structure. Therefore, in the explicit 
description of the iteration procedure, the plane index k 
should be attached to all correlation and auxiliary func- 
tion. 

The logarithm of the cluster functions @' may be in- 
terpreted as effective fields acting on a plane, simulating 
the effect of the half-lattice already summed up. How- 
ever, for simplicity, the functions G' themselves will be 
called effective fields. 

The computational iteration scheme of the cluster 
transfer-matrix method for 3D ANNNI model is as fol- 
lows: 

1. From the cluster functions (effective field) 
k'd'ij^"''''^'' known from the previous step the approxi- 
mate auxiliary function ^k{Sk, Sk+i) = Jlij k&'l^'^'^'^'' 
is constructed and '^k{Sk, Sk+i) in (3) is replaced by it. 

2. The correlation function fe+i$^+"'^+' (12) of the 
auxiliary function ^fk+iiSk+i, Sk+2) is calculated from 
(3). As the both functions at the left-hand side of (3) 
are factorized this problem is equivalent to calculation of 
a correlation function of a 2D lattice model with short- 
range interactions that was discussed above and in previ- 
ous papers 0, 0, 0| in detail. This task is performed 
in two steps, and the approximate factorization utilizing 
ID version of (10) is applied once. 

3. Formula (13) is used and the cluster functions 
fe+iO-'^"'^"*"' are found. 

4. Calculation is continued for the next plane starting 
from the step 1. 

It is convenient to take the result of a previous it- 
eration for the initial condition of the calculation at a 
nearby point in the parameter space. The bulk val- 
ues of the cluster function are obtained after iteration 
over few periods of the commensurate or incommensu- 
rate structure. However, the periods of the commen- 
surate structure sometimes exceed several hundreds of 
lattice constants in our calculations. The convergence 
of the iteration procedure is very slow near the contin- 
uous incommensurate-disorder phase transition and the 
steady state were reached after more than ten thousands 
steps. 

In our actual calculations the length of the cluster 
edges n and I was taken equal to 1, i.e. the cluster on 
which the functions and 6'*^"'"''''' are defined 

has 8 sites (elementary cube) and the functions acquire 
256 values. Thus, our generalized mean-field approxima- 
tion utilizes 256 effective fields instead of one in previous 
approaches. 

The planes perpendicular to the nnn interaction are 
ferromagnetic in the ANNNI model, thus the cluster cor- 
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relation function and the eluster function do not depend 



on its position in the plane. 



is in fact only a 

d+n,j+l 



short-hand notation of ^kikSl^"'"''^ ), where kSl j "'" ' ' is 
a spin configuration of a cluster in the plane k. Similarly, 

To find the actual structure at the given point of the 
phase diagram, it is not necessary to calculate the lat- 
tice site magnetizations. The structure can be deduced 
from the plane dependence of the effective field 0^. In 
our approximation it acquires 256 values, but a plot of 
arbitrary one of them can be used to find the phase di- 
agram. For reason of simplicity and symmetry, the dif- 
ference V'fc = 0fe(+) ^ ©fc(~)' is plotted where "-f" and 
"—"denote spin configurations of the 8-site cluster with 
all the spins up and down, respectively. The sign of the 
function ipk is the same as the sign of the magnetization 
of the A;-th row. The ANNNI model structures consist 
of sequences of planes with negative or positive magne- 
tization. As the external magnetic field is equal to zero, 
the commensurate structures are symmetric with respect 
to spin inversion. Therefore, only is taken into ac- 
count further. Its periodicity p is one half of or equal 
to the structure periodicity if in the interval p the func- 
tion ipk changes its sign even or odd times, respectively. 
A structure consisting repeatedly of p planes with posi- 
tive magnetization and p planes of negative magnetiza- 
tion with periodicity 2p is usually denoted in literature 
as (p). More generally, the sequence of n clusters of the 
above-mentioned planes of the length p interrupted by 
one cluster of the length p — 1 is denoted as {p"{p — 1))- 

At high temperature when the convergence is slow and 
the areas of commensiirate structures arc very narrow, or 
at lower temperatures when near the accumulation points 
the periodicity of commensurate structures tends to in- 
finity, it is often not possible to perform the calculation 
directly at the point of parameter space of desired prop- 
erties, because its precise position is not known. Nev- 
ertheless, the structure at it can be deduced from the 
behavior of the effective field in its close vicinity. For 
this purpose, we shall further plot \ipk+p \ — \ fpk\ vs. \ipk\, 
where p is the periodicity of the function somewhere 
near the point of the parameter space where we perform 
the calculation. It is not necessary to plot |''/'fe+p| — |'0fc| 
vs. \tpk\ for all values of k. The information, we are 
interested in, can be found from behavior of the plot 
for the planes ko + np (n = 1,2, .. .). ko should be the 
number of the plane closest to a node of the structure 
(sign-0fcg ^ sign?/;fe„+i), where [V'/col is close to zero and 
li^kol - IV'feo-il is large. Now, \ipk+p\ - IV'fel =0 always 
means that the structure is commensurate with period p 
or 2p and not that \ipk\ is close to its maximum value. 
The plots will be drawn for in the range from 
to its maximum value when a new plane k — 1 with a 
smaller value of l-fAfe-il appears. Analysis of them will 
make possible to distinguish between commensurate and 
incommensurate structures and confirm the existence of 
the accumulation point, where period of commensurate 
structures tends to infinity. 




FIG. 2: Phase diagram of the 3D ANNNI model. Only the 
basic structures with short periodicity are depicted. The dot- 
ted lines connect points in the parameter space where the 
incommensurate phase has the same periodicity as the corre- 
sponding commensurate structure. The symbols denote the 
following periodically repeating structures; (2) - ttii, (3) - 
TTTiii, (23) - TTiii, (4) - TTTTiiii, F - ferromagnetic (T), 
D - disordered, where the arrows indicate directions of plane 
magnetization. 



III. RESULTS AND DISCUSSION 

Results of our effective-field calculations are consistent 
with the phase diagram obtained by the mean-field ap- 
proximation and low-temperature expansion [sell]. How- 
ever, the temperatures, at which the phase transitions 
occur, are more realistic, and for the exactly soluble case 
J2 = the critical temperature does not deviate more 
than 1% from the true value for the approximation with 
256 effective fields. 

From our calculation, in accordance with previous re- 
sults of other authors, it is possible to conclude that the 
phase diagram consists of infinitely many commensurate 
phases which appear mostly at low temperatures and an 
incommensurate and disordered phase at high tempera- 
tures. 

At low temperature we have found a ferromagnetic 
phase, a commensurate structure with periodicity 4 con- 
sisting of a sequence of couples of planes with alternating 
magnetization ((2)), a structure with periodicity 6 ((3)) 
and combinations of the last two structures of the type 

(2"3) n = 1,2,3, As the low-temperature region is 

fairly well described by the low-temperature expansion, 
we concentrate to the medium- and high-temperature 
properties of the phases (4), (3), (23) and the regions 
in their close vicinity. 

The main phases of the 3D ANNNI model obtained 
from our calculations are shown in the phase diagram 
(Fig. 2). The thick lines denote the borders of the re- 
gions of commensurate phases and represent first-order 
phase transition lines. The dotted lines connect points 
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FIG. 3: Plot of Aijj = \ipk+3\ - IV'fcl vs. tp = [tpkl for every 
third plane of (3) structure. T = 3.45. The plots are drawn for 
the following values of the parameters: 1 - J2/J1 = 0.53335, 2 
- J2/J1 = 0.53325, 3- J2/J1 = 0.53315. The commensurate 
(3) phase is represented by the stable fixed point. Curve 1 
represents an incommensurate structure. As Aip <C ^, the 
plots are practically continuous. 



in the parameter space where the incommensurate phase 
has the same periodicity as the corresponding commensu- 
rate structure. The widths of the commensurate phases 
near the order-disorder phase transition line go to zero 
for all of them, i.e., there is no direct transition between 
the commensurate and the disordered phase. The com- 
mensurate regions at high temperatures are very narrow 
(narrower than the line thickness), nevertheless, they per- 
sist to rather high temperatures. A very large (probably 
infinite) number of commensurate phases between each 
two main phases are not depicted in the diagram and are 
discussed later. The Lifshitz point behind the left edge of 
the diagram is not shown, as the slow convergence of cal- 
culations and complicated phase structure did not make 
possible to correctly interpret the obtained results. 

It is not easy to prove the existence of the commen- 
surate phase in a very narrow region and distinguish be- 
tween the commensurate and incommensurate phase of 
the same or a slightly different periodicity. Here, it is 
helpful to observe the above-mentioned plot of Aip = 
\ipk+p\ — iV'fcl vs. IV'fel, where p is the periodicity of the 
function \ipk\ for the assumed commensurate structure. 

In Fig. 3 this plot for T = 3.45 and two different J2/J1 
inside and near the structure (3) is shown, p = 3, and k 
runs over all planes after which tpk changes its sign, i.e. 
we plot the function for every third plane. The structure 
is commensurate if A'^ is equal to zero. We see that for 
J2/J1 = 0.53335 it never occurs. The function is incom- 
mensurate one with local periodicity greater than 3, and 
it changes at different places of the structure, i.e. the 
true periodicity is very large, and for decreasing J2/J1 
it tends to infinity. It can be considered as a phase- 
modulated (3) structure. As the curve |'0fc+p| ~ IV'fcl vs. 




FIG. 4: Plot of Alp = \ipk+3 \ — \ipk\ vs. V = IV'fcl for every 
third plane and high-temperature, T — 3.52, (3) structure. 
The plots arc drawn for the following values of the parameters: 
1 - J2/J1 = 0.534400, 2 - J2/J1 = 0.534350, 3 - J2/J1 = 
0.534320, 4 - J2/J1 = 0.534313, 5 - J2/J1 = 0.534290. The 
commensurate (3) phase is represented by curve 4. 



I'^fcl shifts in vertical direction with change of J2/J1 with 
only a small change of its shape, we can expect that for 
some values of J2/J1 the curve intersects the x-axis and 
the structure becomes commensurate. In Fig. 3 this sit- 
uation is exemplified by the curves for J2/J1 = 0.53325 
and 0.53315. In the bulk the structure is (3), the differ- 
ence IV'fe-i-pl — \tpk\ is equal to zero and the structure is 
trapped in the stable fixed point. In the transition pe- 
riod, near the surface or a planar defect, where Aip ^ 0, 
the structure is inconimcnsuratc^-like. Starting away from 
an arbitrary boundary condition the system reaches very 
fast an incommensurate metastable state, represented by 
one of the curves, from which the stable bulk commen- 
surate structure at the intersection with x-axis slowly 
develops. 

To confirm the existence of the commensurate struc- 
ture in some region of parameter space, it is not neces- 
sary to find the point where after many iteration steps 
the system converges to a hulk commensurate structure. 
Near to it Aijj is small and the convergence is very slow. 
It is enough to find a nonmonotonous behaviour of the 
Alp vs. ip plot of an incommensurate structure some- 
where near that parameter space point. A set of such 
plots for T = 3.52 is shown in Fig. 4. It is seen that the 
amplitude of modulation of the functions increases when 
approaching the commensurate phase. Thus, already a 
small modulation of the curve far from the commensurate 
structure indicates its presence. 

The width of the region of the commensurate struc- 
ture can be deduced from the amplitude of the function 
I^/i/j+pI — I I vs. \ipk\ and the rate of its vertical shift 
with change of the parameters. 

Fig. 5 shows that for T = 3.546 no commensurate (3) 
structure exists. All the curves are monotonous. The 
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FIG. 5: In close vicinity of the phase transition line to the 
disordered state all the plots A^J = IV'fc+sl ~ \''Pk\ vs. ip = 
\tpk \ are monotonous indicating the absence of commensurate 
phase in this region. 1 T = 3.5461, J2/J1 = 0.5341; 2 - 
T = 3.5465, J2/J1 = 0.5345; 3 - T = 3.5466, J2/J1 = 0.5346. 



■D 




3.56 



FIG. 6: Plot of the width of the commensurate phase (3) vs. 
temperature. The width is given in the units of J2/J1 • 10~*. 



sign of their derivatives is negative and positive above 
and below the x-axis, respectively, and they do not in- 
tersect it. The period of ip close to 3. The point in the 
parameter space is now very close to the order-disorder 
phase transition line, so that the parameters in Fig. 5 
should be carefully changed only in the direction paral- 
lel to it. The rate of convergence is very slow here, and 
the bulk incommensurate structures depicted in the fig- 
ure were obtained after more than 10,000 iteration steps. 
The argument of nonexistence of commensurate struc- 
ture at T = 3.546 is confirmed by extrapolation of the 
widths of the (3) structure to higher temperatures shown 
in Fig. 6. For small values of the width, this plot could 
be well fitted by a parabola. 

Similar considerations were done for the structures (4) , 
(23) and (23^") and it was found that the commensu- 
rate structures of higher periodicity disappear at lower 
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FIG. 7: Number n of ttt or HI piano sequences in (23") 
commensurate phases near the transition line to (3) structure. 
The value of J2c at the transition line is different for each 
temperature, n acquires discrete values and the lines are only 
guides for the eye. 



temperatures. The whole region near the order-disorder 
phase transition line is incommensurate with tongues of 
commensurate structures of low periodicity which do not 
reach the phase transition line. 

At very low temperatures the phase (3) is neighboring 
to the phase (23). With increasing temperature, at T = 
1.6, phases of the type (23") start to appear. At given 
temperature T the period of the structure, p = 3n + 2, 
increases with decreasing J2/J1 and its largest value is 
reached at the boundary of the (3) structure. The plots of 
n near the (3) boundary for T = 2.60, 2.62, 2.63, 2.64, 2.65 
are depicted in Fig. 7. The periodicity in the close vicin- 
ity of (3) phase increases very fast and for the tempera- 
tures above 2.62 it is not possible to determine the value 
of Timax or even decide if it is finite or not. Nevertheless, 
the accumulation point, where rimax becomes infinite, can 
be found analyzing the plots lipk+p] — I^Afel vs. I'i/'fel for dif- 
ferent temperatures and p = 3, which are shown in Fig. 8. 
The structure (23") for large n is formed from domains 
of the structure (3) of the length slightly less than n in- 
terrupted by domain walls symbolically denoted by "2" 
in the symbol (23"). The (3) structure beyond the wall 
is shifted by one plane with respect to the structure in 
the previous domain. The (3) domains correspond to the 
minima of the plots in Fig. 8 where the iV'ft-i-sl — iV'ftI are 
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FIG. 8: Plot of Alp = lipk+s] — IV"*;! vs. ip = \tpk\ for com- 
mensurate (23") structures near transition to (3) phase, n is 
equal to the number of points along each curve. (Only points 
for small value of Aip are depicted in the figure.) The curves 
for J2 = J2c touch the x-axis. For curves with zero derivative 
in the minimum n tends to infinity. 



practically equal to zero. The advent of the wall is so 
abrupt and IV'fe+sl — \4'k\ so large that the next point af- 
ter the very right edge of the each curve is already out of 
scope of the diagram. With decreasing J2/J1 the plots 
are shifting down and at the (3) phase transition line the 
minimum of the plot touches x-axis, and after some tran- 
sition period the system remains stuck in (3) phase. If 
the slope of the plot at minimum is zero, the periodicity 
near the boundary tends to infinity, and the temperature 
of the system is already above the accumulation point. 
From Fig. 8 we see that the accumulation point is close to 
the temperature T — 2.64. The periodicity of tends 
to infinity if the curve approaches x-axis for J2 —>■ J2C 
a.t T — 2.65. Using our method, we were able to find a 
commensurate structure ofp = 1802 at this temperature. 

For temperatures lower than the temperature of accu- 
mulation point near the boundary of the (3) phase only 
plain (23") phases exist. This is in contradiction with 
the simple mean field approximation findings 01 ■ There 
is also a difference in location of the accumulation point. 
In 1^ it was found well below the turning point of the (3) 
boundary where the width of the (3) phase starts becom- 
ing narrower. Our approach locates the accumulation 
point T — 2.64 slightly above the turning point. Approx- 
imately at the same temperature first combined phases of 



the type ((23")(23"+i))) between ((23")) and ((23"+i))) 
for large n start to appear. Similarly to the previous 
more simple case, following combinations in the hierarchy 
arc of the type (((23"))'''(23"+i))) or ((23")((23"+i))'=)), 
which appear at temperatures by 0.04 higher than the 
temperature of the accumulation point. A great compu- 
tational effort is needed to detect a next type of combina- 
tions ((((23"))'=(23"+i)))'((((23"))*=+i(23"+i)))}. These 
high-order combinations occupy very small areas of the 
parameter space, and with increasing temperature they 
are soon replaced by incommensurate structures. 

It is widely believed that the structures with large 
distances between the domain walls are commensurate 
whereas the structures where the distance between the 
walls is shorter than the wall-wall interaction are incom- 
mensurate. This statement should be formulated more 
precisely. The commensurate structures with short dis- 
tances between the walls are more stable than those 
with longer distances. They persist to higher temper- 
atures and they occupy a wider area in the parameter 
space. Nevertheless, in the areas between ((23")) and 
((23""''^))) for small n the onset of incommensurate struc- 
tures was found at lower temperatures than for large n. 
The distance between the commensurate structures in- 
creases with decreasing n faster than their width so that 
there is enough space for incommensurate structures. 

In summary, we developed an effective field approxima- 
tion, wihich yields by simple iteration procedure practi- 
cally any of probably an infinite number of phases in the 
phase diagram of the 3D ANNNI model. In fact, the 
method treats an infinite lattice. Lattice size, in con- 
trast to other mean-field and DMRG approaches, does 
not enter the calculation. A difficult task to distinguish 
between commensurate and incommensurate structure 
after a finite number of iteration was made easier by 
plotting the derivative of an effective field with respect 
to its value in course if iteration. Our calculations con- 
firmed the general picture of the phase diagram obtained 
by other methods, made it more accurate and supported 
the suggestion following from the Monte Carlo calcula- 
tions that the commensurate phases are separated from 
the disordered phase by an incommensurate region. 
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